//try a little hourly wage regression
use "$temp/psid_long_cleaned", replace
gen exp = age-(educ + 6)
gen wage_hourly = wages/hours
gen lwage = log(wage_hourly)
//gen coll = (educ>=16)
keep if !coll
areg lwage c.exp##c.exp##c.exp##c.exp educ ib19.statefips i.year [aw=weight] , cl(uniqid) a(uniqid)


//convert to data
levelsof statefips, clean local (fips)
clear
set obs 51
gen statefips = .
gen skill_price_psid = .
gen se = .
local counter = 0

foreach fip in `fips'{
	local counter `++counter'
	
	replace statefips = `fip' in `counter'
	replace skill_price_psid = _b[`fip'.statefips] in `counter'
	replace se = _se[`fip'.statefips] in `counter'
}	

//exponentiate and save
gen ub = skill_price_psid + 1.96 * se
gen lb = skill_price_psid - 1.96 * se


replace skill_price_psid = exp(skill_price_psid)
replace ub = exp(ub)
replace lb = exp(lb)
drop se

save "$temp/skill_prices_psid", replace
merge 1:1 statefips using "$temp/state_dynamics_all", keep(match) nogen
corr skill_price_psid skill_price_2000
scatter skill_price_psid skill_price_2000 || line skill_price_2000 skill_price_2000, bgcolor(white) graphregion(color(white)) ytitle("PSID Skill Price") xtitle("Baseline Skill Price") legend(off)
graph export "$output/skill_price_robustness_psid.png", replace
graph close

count if skill_price_2000 < lb | skill_price_2000>ub




//